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Introduction. The paper deals with the calculation of wooden arches taking into account the nonlinear relationship be- 
tween stresses and instantaneous deformations, as well as creep and geometric nonlinearity, are considered. The analy- 
sis 1s based on the integral equation of the viscoelastoplastic hereditary aging model, originally proposed by 
A.G. Tamrazyan [1] to describe the nonlinear creep of concrete. 

Materials and Methods. The creep measure is taken in accordance with the work of I.E. Prokopovich and 
V.A. Zedgenidze [2] as a sum of exponential functions. The transition from the integral form of the creep law to the 
differential form is shown. The relationship between stresses and instantaneous deformations for wood under compres- 
sion is determined from the Gerstner formula, and elastic work is assumed under tension. The solution is carried out 
using the finite element method in combination with the Newton-Raphson method and the Euler method according to 
the scheme that involves a stepwise increase 1n the load with correction of the stiffness matrix taking into account the 
change in the coordinates of the nodes with the sequential calculation of additional displacements of the nodes, which 
are due to the residual forces. The proposed approach for increasing the accuracy of determination of creep defor- 
mations at each step provides using the fourth-order Runge-Kutta method instead of the Euler method. 

Results. Based on the Lagrange variational principle, expressions are obtained for the stiffness matrix and the vector of 
additional dummy loads due to creep. The method developed by the authors is implemented in the form of a program in 
the MATLAB environment. Calculation examples are given for parabolic arches simply supported at the ends without 
an intermediate hinge and with an intermediate hinge in the middle of the span under the action of a uniformly distrib- 
uted load. The results obtained are compared in the viscoelastic and viscoelastic formulation. The reliability of the re- 
sults is validated through the calculation in the elastic formulation in the ANSYS software package. 

Discussion and Conclusions. For the arches considered, it 1s found that even with a load close to the instant critical, the 
growth of time travel is limited. Thus, the nature of their work under creep conditions differs drastically from the nature 


of the deformation of compressed rods. 
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Introduction. Wood refers to the materials that exhibit their non-linear properties, both under short-term and 
long-term exposure. Most existing rheological models of wood establish its instantaneous properties on the basis of 
Hooke's law; but for compressed wood, a nonlinear diagram of the elastoplastic type is typical [3]. 
K. P. Pyatikrestovsky was the first to investigate the issues of joint accounting of instantaneous elastoplastic properties 
of wood and its creep. He combined these properties on the basis of the method of modulus of long-term deformation 


[4-7]. In this method, the creep equation contains time explicitly, which significantly limits its use in the case of 
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complex loading modes and variable loads. There is a significant number of publications on the joint accounting of 
creep and instantaneous nonlinearity of wood deformation in the calculation of compressed rods [8-14]. In particular, 
[8] discusses the applicability of the Boltzmann superposition principle to the description of nonlinear wood creep. In 
[9], theoretical calculations are compared to experimental data for compressed rods; and in [10], in addition to the 
compressive longitudinal load, the transverse load is taken into account. In [11—13], the stability problem is solved for 
the rod elements in the frame structure. In paper [14], a long-term critical force is derived for compressed wooden 
elements, taking into account the nonlinear creep. 

The problem of calculation with account for creep and instantaneous nonlinearity of deformation 1s a challenge 
not only for any particular rod and frame, but also for such rod systems as arches. This work objective is to develop a 
method for calculating arch structures taking into account the nonlinear properties of the material under short-term and 
long-term exposure, as well as geometric nonlinearity. 

Materials and Methods. As a relation that determines the relationship between stresses and 


deformations, we use the equation of the viscoelastic model of hereditary aging: 


6C(t,T) 


(1) 
(A. jrfo(ry] ear, 


r 


This equation was first proposed in [1] for modeling the nonlinear creep of concrete. Here, the function f(o) 
establishes the relationship between stress and instantaneous strain, C(¢,r) — creep measure. The stress-strain diagram 


of compressed wood under short-term loading is well approximated by the Gerstner formula [15], which has the form: 


2 


E 
o=Fe-7 oe. (2) 


The compressive stresses are used in the formula (2) with the sign “+”. When wood 1s stretched, there is a 


linear diagram up to the breaking point. Expressing the strain in terms of stress from (2), we obtain: 
es | i saa | 
E, R 
Based on (3), the stress function f(o) can be written as: 


f(o)= 2e{i-i-£ 
(4) 


To measure the creep of wood, we use the formula proposed in the work of V. A. Zedgenidze and 
I. E. Prokopovich [2]: 


(3) 


C(t,t) = ce es Ae" )pi-Be™™ J, (5) 
where Cy = 2.87:10° MPa”, Ay = 10.95-10° MPa’, B, = 1, y=, =0.15 day 7. 
Equation (1) can be presented in the form: 


Ed 


Oo 
S=—+E., 6 
E (6) 


dt — creep deformation. 





Eo 6C (t,t) 
here E = —°— — secant modulus, ¢* =— pat a 2 
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For the creep measure in the form (5), the creep strain can be written as the sum of two components: 








6 =e +6,, & = 1C,B, | f.o(z)] e“ dq, = Ay |f{.o(z)] aU: (7) 
Differentiating in time (7), we obtain expressions for the growth rates of each component: 
ae’ Oe . 
= =7(CBF[o(t)]-e}); == Arf Lo(s)]e”. (8) 


On the basis of (6), we derive the relationship between internal forces and deformations for the element under 
the combined action of the longitudinal force and the bending moment. The elastic modulus is assumed to be a function 
of y coordinate, which varies from -//2 to h/2 in the height of the cross section. Based on the hypothesis of flat sections, 


we write the total deformation in the form: 


Mechanics 


115 


Advanced Engineering Research 2021. V. 21, no. 2. P. 114-122. ISSN 2687-1653 





E=6,+ VX, (9) 


2 





Vv 
— change in curvature.. 


where &) — deformation of the middle layer ¥ =—— 
bs 


We then substitute (9) in (6). Expressing o by €, we obtain: 
o=E(y)(e-e' )=E(y)(e,+yx-€ ). (10) 
The longitudinal force and bending moment in the element are related to the stress through the following 


integral dependences: 


M =|oydd; N= | odd. (11) 
A A 


Here, A — cross-sectional area of the rod. 


Substitute (10) in (11) and convert the resulting equalities to the matrix form: 


{aa} FL (12) 


: : , : EA ES ; 
where N =[E(y)e dA, M =[E(y)ye dA , PI=| es a — the matrix of reduced stiffnesses, which are 
determined from the formulas: 

EA= JE») ) dA, ES = JE(y ) ydA, EI = JE(y) ) dA. (13) 


The problem with account for the physical and geometric nonlinearity will be solved by the authors using the 
finite element method. The rod element shown in Fig. 1 is used. The axial strain, taking into account the geometric 
nonlinearity, is the sum of the linear and nonlinear components: 

du 1_dyv,., 


&, =6, +6 =—+=— 14 
0 0+ ene (14) 





XxX: 


~ 


J 
Fig. 1. Rod finite element 


To obtain a system of FEM equations, we use the Lagrange variational principle. The strain potential energy 
(SPE) is written as: 
1 ; 1 : my 
I == | oe"dV = 5 E(v)(e +8 + K-82 ) dV. (15) 
V V 


where ¢° — elastic strain, equal to the difference between total strain and creep strain. 


We write (15) as the sum of four integrals: 

n=1{ Jee (E, +60 + yy-€ av + |B(y) ef (e, +e, + y—€ )dV+ 
(16) 
+] E(y) vx(e +6" + yy—e')dV-[E( ye (, +6 ns )aV| 
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The first integral in (16) is written as: 
I 
JE(y) (E, +8" +yy-€ )dV = |e, [E(y)(e), + yx-€' )dddx+ (17) 
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“fe 81 | E(y) dAdx = Jena fo .¢". EAdx, 


where / — the length of the finite element. 


The second integral in (16) is presented in the following form: 


JE(y) (E, +E + yy-€ )dV = fe € J ECy ) dA | EC y)( (e") dV+ 
(18) 
+] e} -x[E(v) ydddx— |e; | B( ye’ dAdx, 
The second term in the right-hand side of (18) can be neglected, given its higher order of 
smallness. The third integral in (16) is written as: 
JE(») yx é, te, + yy-€ )dV = jxfAlv ) We, + yx, —€' ) dAdx+ (19) 
l l l 
+ fe, ESdx = [mu yd + [ery - BSdx. 
0 0 0 
The fourth integral in (16) is presented in the following form: 
JE(») (E, +6 + yy-€ )dV = ae ) ye dAdx+ 
l l 
+fe, [e'E(y)d4dx—|(e') E(y)aV + #5 fe E(y) ) dAdx = (20) 
0 A V 0 


l l l 
= [yM dx + [ede + [et N’dx-[(e°) E(y)aV 
0 0 0 V 
The summand | (ey E(y)dV when minimized over the vector of nodal displacements, will vanish. Finally, 


the expression for the SPE will take the form: 


i I I i I l 
i= = é, Ndx + 2 6, EAdx + | Mya + 2| ei yESdx —|y¥M" dx -| e, N dx + 
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l . 1 l N N- l 
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where N, =¢,FA+y,,ES—N° — the average axial force in the element, y,, — the mean change in the curvature of 


(21) 


the element. 


Taking into account (12), formula (21) will take the form: 
1 l N- 
T=—| fs" | [D]{e}-2 dx+N _ |e dx , 22 
+[ te’ [[evte)-2{ bacon, fe an 


For finite element displacements, we assume the following approximation: 


where {e}’ ={e x}. 


u, —u 





u(x) =u, + _ x; (23) 
v(x) =O, +0,xX+0,X° +0,x°; (24) 
p(x) = - =O, +20,x+3a0,x°. (25) 


The coefficients of polynomial (24) are determined through substituting the coordinates of the following nodes 
in (24) and (25): 
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From (26), vector {a}={a, a, ©, @,}’ is expressed as follows: 
0 0 0 0 U, 
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0 1 ®, . 2 3 1}) % 
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Taking into account (27), formula (24) will take the form: 
={V}{U}. (28) 
roe M={1 x x wh]. 
Vector {¢} is written as: 
= -1/I 00 1/1 00 
x ; = 
t=] [=| a Ph Uy = [By (29) 
ie dx 
After substituting (29) and (28) in (22), the SPE is written as: 
dip’ dy¥ 
=~ Uy ill D||B]dx{U uy flay “\M peor N Z 3 — 7 ae(U}. (30) 


After minimizing the Lagrange functional with respect to the nodal displacement vector, we arrive at a system 


of equations of the following form: 
+{Fh, (31) 
APY atv 

dx 


([K]+[K, |){U} ={F} 


I I 
where [K ]=|l47 |D|[B]dx — stiffness matrix, | K, |=N,, | ba — geometric stiffness matrix, 








l N° 
= [t BY dx te — contribution to the load vector of creep strains, {F'} — vector of external nodal forces. 


A physically and geometrically nonlinear problem is solved using the Newton-Raphson method. The first step 
is the calculation at t = 0. The load increment is performed in quasi-statically small portions with the sequential 
calculation of additional displacements of the nodes, which are caused by the residual forces, the adjustment of the 


tangent modulus of elasticity and the coordinates of the nodes at each step. 








Research Results. The presented equations and the calculation algorithm are implemented by the authors in 


~ 

= Next, the time interval at which the calculation is performed is divided into a finite number of time steps. The 
— 

: creep calculation is performed in the same way as the static load calculation. The increment of creep deformations at the 
> step ttAt can be determined using the Euler method 

~~ 

(s . (es maz 

7; Ag = > | At. (32) 
o ar Ot 

= 

ir. 

= 

on 


the form of a program in the MATLAB environment. The calculation of a parabolic arch simply supported at the ends 
under the action of a load evenly distributed along the length was performed (Fig. 2) with the following initial data: Ep = 
1.48-10° MPa, R = 55 MPa, L = 16m, f=3.2 m. The cross-section of the arch was assumed to be rectangular with 


118 dimensions 6 = 10 cm, 4 = 15 cm. The arch was divided into 40 finite elements in length, the section was divided into 


Yazyev S. B., et al. Stability analysis of wooden arches with account for nonlinear creep 





100 segments in height, the number of steps in time was assumed to be 600, and in load — 200. The maximum number 
of iterations at each step is 20. Figure 3 shows the maximum deflection versus load graph. The dashed line corresponds 
to the calculation in the elastic formulation. To check the results, the solution to the elastic problem was also performed 
in the ANSYS software package using BEAM 188 rod finite elements. The number of finite elements was assumed to 
be the same as in MATLAB. There was no significant difference in the results. When calculating in the elastic 
formulation, a sharp increase in displacements corresponding to the loss of stability 1s observed at g = 10.5 kN/m, and 


when taking into account the instantaneous nonlinearity of deformation — at g = 10 kN/m. 


q 


SERRA 


L 


Fig. 2. Design scheme of the structure 


gq — distributed load on the arch, f— arch rise, L — arch span 
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Fig. 3. Dependence of maximum arch deflection on the load under short-term loading 


The creep calculation for this arch has shown that even at loads close enough to the instantaneous critical load, 
the growth of displacements is limited. Figure 4 shows a graph of the midspan deflection growth versus time under a 
load of g = 8 kN/m. The dashed line on this graph corresponds to the solution in the viscoelastic formulation. There was 
no significant difference in the results. 

Also, with the above initial data, the calculation of a three-pinned arch with an intermediate hinge in the 
middle of the span was performed. In this case, the instantaneous critical load was significantly lower. When calculated 
in the elastic formulation, it was 4 kN/m, and with account for the instantaneous nonlinearity of deformation — 
3.3 kN/m. As in the previous example, even when the load is close enough to the instantaneous critical, the creep 
decays. The curves of change in time of maximum deflection at g = 3 kN/m are shown in Fig. 5. The symbols are the 
same as in the previous graph. 


Mechanics 


119 


http://vestnik-donstu.ru 


120 


Advanced Engineering Research 2021. V. 21, no. 2. P. 114-122. ISSN 2687-1653 
SSS —————————— a 
w, mm 


24 


22 
20 
18 
16 
14 
12 
10 
8] 


0 50 100 150 200 250 300 
t, day 


Fig. 4. Change in deflection in the middle of the span over time under load g = 8 kN/m 


—_ 
—_— = 





0 50 100 150 200 250 300 
t, day 


Fig. 5. The growth in time of maximum deflection for an arch with an intermediate hinge in the middle of the span at q = 3 kN/m 


Discussion and Conclusions. The resulting equations and the developed method are universal and provide the 
use of arbitrary dependences between stresses and instantaneous strains, as well as arbitrary expressions for the creep 
measure. This provides calculations of structures made not only of wood, but also of any other material. As a result of 
the analysis of the creep of wooden arches, it is found that, in contrast to compressed rods, creep is limited for them, 


even under loads close to the instantaneous critical load. 
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